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We present a model intended for rapid sampling of ground and excited state potential energy 
surfaces for first-row transition metal active sites. The method is computationally inexpensive and 
is suited for dynamics simulations where (1) adiabatic states are required "on-the-fly" and (2) the 
primary source of the electronic coupling between the diabatic states is the perturbative spin-orbit 
interaction among the 3d electrons. The model Hamiltonian we develop is a variant of the Anderson 
impurity model and achieves efficiency through a physically motivated basis set reduction based on 
the large value of the d-d Coulomb interaction Ud and a Ldnczos matrix diagonalization routine 
to solve for eigenvalues. The model parameters are constrained by fits to the partial density of 
states (PDOS) obtained from ab initio density functional theory calculations. For a particular 
application of our model we focus on electron-transfer occuring between cobalt ions solvated by 
ammonium, incorporating configuration interaction between multiplet states for both metal ions. 
We demonstrate the capability of the method to efficiently calculate adiabatic potential energy 
surfaces and the electronic coupling factor we have calculated compares well to previous calculations 
and experiment. 



I. INTRODUCTION 



First-row transition metals (TM) complexed with or- 
ganic ligands play key roles in biology and aqueous chem- 
istry. Ligand exchange and electron transfer are po- 
tentially rate-limiting processes that must be consid- 
ered in industrial or bio-remedial chemistry |l|. Ligand 
binding at the iron porphyrin active sites of heme pro- 
teins are fundamental for oxygen transport and aerobic 
metabolism 0. Many of these reactions involve non- 
radiative transitions between different spin states and 
require spin-orbit coupling. This interaction occurs pri- 
marily on the TM atom due to the localized nature of the 
3d valence electrons and most quantum chemistry codes 
do not incorporate the small (^0.05-0.1 eV) energy-scale 
coupling. 

In general, treatment of open-shell TM atoms is dif- 
ficult and many of the standard density-functional the- 
ory (DFT) model chemistries such as the local density 
or gradient-corrected approximations (LDA and GGA, 
respectively) arc known to give erroneous qualitative re- 
sults in some cases |3| . The primary source of these prob- 
lems are non- dynamical correlations associated with the 
near-degeneracy of configurations with the Hartree-Fock 
single-Slater determinant groundstate |J]. In fact, it is 
clear that in the vicinity of the transition state of a charge 
transfer complex such non- dynamical correlations may 
be very important, much as they are in chemical bond 
breaking. One may incorporate these effects by higher- 
order corrections such as Moller-Plesset perturbation the- 
ory (MP2,MP4) or configuration interaction (CI) meth- 
ods such as complete active space self-consistent field 
(CASSCF), or, especially for non-dynamical correlations, 
multireference configuration interaction approaches [5|, 



but these methods are computationally prohibitive for 
TM molecules containing ^50 atoms. In many cases the 
CI must be limited to single excitations, which may be 
insufficient [6|. 

While the addition of some amount of non-local ex- 
change used in the hybrid Hartree-Fock/DFT methods 
such as B3LYP seems to capture some of the multi- 
determinantal character of TM systems 0,|al3' the min- 
imal local Gaussian basis sets required for quantitative 
accuracy are very large and the sampling required for 
dynamics and the calculation of excited states is not yet 
affordable for most processors. 

In this paper we introduce an effective semi-empirical 
model Hamiltonian that is well-suited for (^50-100 
atoms) TM molecules. The problematic aspect of 
strongly correlated TM systems, the (d-d Coulomb inter- 
action) , leads to a physically- motivated reduction of the 
many-body basis set required by limiting the number of 
d n configurations that must be considered. We develop 
a variant of the Anderson Impurity Model (AIM) , which 
has had some success in describing spectra and energet- 
ics of TM solid-state systems Qij . We then describe our 
method of solving for the groundstate of this model us- 
ing the variational ansatz for the groundstate wavefunc- 
tion of Gunnarsson and Schonhammer j 1 11 ] and Varma 
and Yafet [12] . The AIM reduces the multiple-degrees 
of freedom of this problem to a small number of one- 
electron and Coulomb input parameters. We show how 
these parameters can be constrained and, in some cases, 
directly extracted from ab initio DFT calculations. We 
have done this using the SIESTA yjj code. The tensor 
product basis formed from the TM 3d and ligand orbitals 
results in a very large eigenvalue problem. We solve for 
ground and excited states with an efficient Ldnczos diag- 



onalization routine. This model, which is similar to the 
angular-overlap model of Gerloch and co-workers |14J . 
differs in both the physically-motivated truncation of the 
basis set and also in the method of parameterization. We 
seek here to identify the microscopic origin of model pa- 
rameters that describe the electronic structure, such as 
the d-d Coulomb repulsion and the ligand-field splitting, 
rather than focusing on structural sensitivities 15] 

We now develop the model in detail and apply it to the 
electron transfer reaction between Co(NH 3 ) 2 f , [Co(II)] 

and Co(NH 3 )^ ' , [Co(III)] molecules in solution, which 
has been extensively investigated by ab initio electronic 
structure methods [l|| UA 03 ■ 



II. THE MODEL 



The TM molecules we wish to investigate are typi- 
cally composed of a single 3d TM atom bound to one 
or more ligands. We express the model Hamiltonian in 
2nd-quantized notation using the 3d atomic orbitals of 
the TM atom and the molecular orbital eigenstates of 
the ligands. We have recently applied a similar model 
to the study of the electronic structure of cobalt valence 
tautomers [l9J . 

The electron correlation effects that we wish to incor- 
porate into our model result from the localized character 
of the atomic 3<i-orbitals on the TM atom. The ligand 
molecular orbitals will tend to be more delocalized over 
several atomic sites. A set of localized atomic levels hy- 
bridized to a reservoir of more extended electronic states 
is well-described by an Anderson Impurity Model. We 
may write down a variant of this model, specific for TM 
molecules: 



H = e L Y^ c f ja c ia + Y^^dja - g 1 x]dl fa d ia (1) 
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where en are the on-site energies for the redox-active 
ligand molecular orbitals and c- are the corresponding 
fcrmion creation operators for the ligand molecular or- 
bitals. The operator dl a creates an electron in the TM 
3d-orbital with Oh symmetry label 7 = xy, yz, xz, 3z 2 — 
r 2 ,x 2 — y 2 and a spin. V 1 i are the hopping matrix ele- 
ment between the Co 3<i-orbitals and the ligand states. 
We also include the 3d spin-orbit (s.o.) in its first-order 



perturbative form. The inclusion of the TM s.o. is signif- 
icant for processes that involve changes in the total spin 
of the electrons in the 3d shell. 

Given the large changes in reorganization energy that 
can accompany these processes, we need to couple our 
model, which deals only with the valence electronic sec- 
tor, to a relevant reaction coordinate, such as bond length 
change. This is especially vital for biomolecules where 
function is modulated by structural changes |20|. The 
simplest way to couple the TM molecule to the environ- 
ment is via the terms in Eqn.(4): The term g 1 xdl c ,d lcr 
provides a (Holstein) linear coupling of the d-level charge 
occupancy of the TM atom and the reaction coordinate, 
the change in metal-ligand bond distance, x. At low- 
est order, the potential energy surfaces are harmonic and 
so we can add an energy term that is parabolic in x: 
(l/2)Kx 2 . The linear coupling and the parabolic energy 
term can both be obtained from a gradient expansion of 
the hybridization about the M-L bond length minimum 

m. 
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where V(x) \ x=xo = V ia , the linear coupling -^ lx=xg 
we estimate by taking the first derivative of the ap- 
propriate Slater-Koster i nteg ral |22fl . using a distance- 
parameterized expression |23|. 
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where xq is the equilibrium bond length. This expres- 
sion should be valid for small displacements away from 
xq and weak coupling between the TM d-electrons and 
changes in M-L bond length. From a simple harmonic 
oscillator model K is related to the single-active metal- 
ligand stretching mode ui by u = y/K/M where M is the 
reduced mass which may be approximated by the total 
mass of all six NH3 ligands. To obtain K we refer to 
experimental metal-ligand vibrational data |24| . 

We now show how we apply this model to a particular 
system, the self-exchange electron-transfer in the cobalt 
hexaammines, Co{NHz)J~ a ' . . 



III. ELECTRON-TRANSFER IN Co{NH 3 f+^ ) + 



A. Experimental and Theoretical Review 

We would like to examine how effects arising from 
strong correlations manifest in the general problem of 
electron transfer between solvated l s *-row transition 
metal ions that are mediated by a ligand bridge. In par- 
ticular, we focus on the electron self-exchange rate in 



the cobalt hexaammine system Co(NH 3 ) 6 tJ ■, , 
given by the following reaction: 



Co(NH 3 ) 2 + 



Co(NH 3 ) 3 ^ 



which is 



(4) 



Co(II) 

2 E 



Co(III) 
3 T, 



6(05) ^ ""V" -"-iJeiaq) ,_ 

Co(NH 3 ) 3 6 -t aq) +Co(NH 3 )l+ aq) 

and the reasons why k et is two orders of magnitude larger 
than predicted. Experimental measurements of the rate 
(k et ~10 _6 -10- 5 M -1 s -1 ) [H deem the reaction adi- 
abatic, i.e. where the electron transmission coefficient 
k w 1, while previous theoretical efforts (k et ^f0~ 8 
M _1 s _1 ) yield n e i « 1 J2{|. Electronic structure effects 
can be important for ET between small TM complexes in 
solution (which we define to be the metal ion and the six 
coordinating ligands of the first solvation shell) that are 
in contact with each other. In the case of the hexaam- 
mines, the salient feature ap pea rs to be a "spin barrier", 
(first pointed out by Orgel 27]). During the course of 
the transfer process, each cobalt ion must change its spin 
by AS = 3/2. Strictly speaking, k et should be zero, but 
various schemes have been proposed to show why the rate 
is, in fact, almost adiabatic. 

The model for the electron-transfer process is usually 
treated as two cobalt ions (one Co(II) valence, the other 
Co (III)) that are each surrounded by their primary solva- 
tion shells which consist of six octahedrally-coordinated 
NH 3 groups. The ET between the Co(NH 3 )g +/3+ com- 
plexes is generally believed to be an outer-sphere pro- 
cess, so the ligands act as a bridge that the tunnelling 
electron delocalizes over during the reaction. The redox- 
active orbitals of the ammine ligands are the highest 
occupied molecular orbitals (HOMOs) of the ammincs 
which have a symmetry and hybridize strongly with the 
cobalt e g orbitals |23. The strong a-donating ability of 
the ammine ligand has been confirmed by NMR stud- 
ies which observe the kinetic stability of this complex 
through the rate of water exchange of Co(NH3) 6 
This rate for Co(NH 3 )6 +/3+ is a factor of -10 6 slower 

than for Co(H 2 0)6 +/3+ Q. 

The ET reaction can be reduced to a unimolecular pro- 
cess by consideration of the encounter complex, which 
is formed by the van der Waals contact of the two pri- 
mary coordination shells of the metal ions. This can 
be done in a variety of orientations 29]. The two most 
commonly considered are the (1) apex-to-apex, which in- 
volves contact along the C4 axis, and is the most opti- 
mal bridge since the Is orbitals of the hydrogen atoms 
of the two ammines in close-contact show large mixing, 
and (2) face-to-face, which has two faces perpendicular 
to the C 3 body-diagonal. The ammine contact is proba- 
bly not as good here but there are 6 ammines in contact 
instead of 2. For calculation of k et and Hjja, some aver- 
aging over different orientations should be done to give a 
more realistic estimates of these quantities. We do not do 
that here, only focusing on the most optimal geometry, 
apex-to-apex. Given this fixed geometry, we investigate 
electronic structure effects. 
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FIG. 1: Lowest energy configurations and excited states for 
Co 2+ and Co 3+ valences. These are the relevant electronic 
states for the self-exchange electron transfer reaction between 
cobalt atoms in Co(NH3) 6 . The question mark for the 

Co 2+ excited state refers to an estimate by Buhks, et al. at 
the equlibrium configuration. 



Buhks et al. [2|| first noted that the groundstate- 
to-groundstate electron-transfer between the low-spin 
Co(III) ion (3d 6 )( 1 Ai) and the high-spin Co (II) ion 
(3d 7 )( 4 Ti) is spin-forbidden and the reaction must pro- 
ceed through the S = 1/2 2 E excited state for Co(II) or 
the S = 1 3 Ti and 3 T 2 of Co(III). They also showed that 
there is a non-zero s.o. coupling matrix element between 
these ground and excited states, so the groundstates 
are not single-determinant but show small admixture to 
these excited states. Thermal population of these excited 
configurations is an open issue [la. 1171 fig] . The vertical 
energy separation between x Ai and 3 Ti for Co(III) is 
very large (~14,000 cm -1 ) and well-established by ex- 
periment. On the other hand, the separation between 
4 Ti and 2 E for Co(II) has been estimated from semi- 
empirical INDO/CI calculations to be 3000-9000 cm" 1 
|16| . While the existence of low-spin Co(III) has been 
experimentally verified, no such confirmation exists for 
high-spin Co(II). Thus, appreciable population of the 2 E 
may yet be a viable possibility. Also, it should be noted 
that these values have been estimated for the equilibrium 
configuration. The relevant values for the energy separa- 
tion should be made at the transition state, where they 
will be significantly reduced. This discussion is summa- 
rized in Fig.(l). 

Other theoretical work, most notably by Newton |2S| 
and Sutin 30], who used ab initio self-consistent field, 
spin-restricted Hartree-Fock to calculate the electronic 
coupling matrix element Hda, assumed 2 E as the Co(II) 
groundstate. Newton obtained a value of 940 cm -1 for 
the spin-allowed electron transfer. He then used a reduc- 
tion factor of 2 X 10 -2 to account for the 2 nd -order spin- 
orbit coupling, obtaining a value of 20 cm -1 for the 4 Ti 
to : Ai process which corresponds toa«~ 4.0 X 10~ 3 . 



Furthermore, Newton considered different orientations 
for the encounter complex, such as apex-to-apex, edge- 
to-edge, and face-to-face, which will alter hybridization- 
matrix elements and differ in the number of bridging lig- 
ands (2, 4, 6, respectively). 

Newton, in a more recent paper [18|, calculated the 
relative energies of 4 Ti and 2 E electronic configura- 
tions of Co(II) and the 1 A i , 3 Ti of Co(III) using un- 
restricted Hartree-Fock (UHF) and second-order M0ller- 
Plesset (UMP2). He recalculated Hda and k for the 
apex-to-apex geometry using these semi-empirical CI 
methods that take into account electron correlations at 
the level of 2nd-order Brillouin-Wigner perturbation the- 
ory. For the groundstate pathway 4 Ti/ 1 Ai ^ 1 Ai/ 4 Ti 
he obtains a values (Hda — 9.5 cm -1 , n — 1.2 X 10~ 3 ) 
which are optimal in the sense that the 3z 2 — r 2 - NH3 
(T-HOMO hybridization can only decrease in other ge- 
ometries. He estimates the effect of orientational aver- 
aging by taking into account a multiplicative factor 0.1; 
K e i ~Q.l(n e i)apex-to-apex so K e i ~10~ 4 , still falling 2 or- 
ders of magnitude short of explaining the observed re- 
sults. 



B. Many-Electron Basis Set 

We use the model Hamiltonian of Eqn.(l), which in- 
volves use of an ansatz for the groundstate wavefunc- 
tion with the coefficients to be determined variation- 
ally. The wavefunction is composed of the lowest-lying 
single-determinants in the many-electron Hilbert space, 
where the truncation of the space is dictated by energet- 
ics. Our problem therefore reduces to diagonalizing the 
model Hamiltonian in a restricted Hilbert space. We first 
considered the two nominal valences for each cobalt hex- 
ammine reactant. In each case the most significant states 
are the high-spin and low-spin state, the spin-orbit ex- 
cited states, and states that couple by hole transport to 
the bridging NH 3 HOMO. For Co(NH 3 ) 2+ we use the 
groundstate 4 Ti, which is the three-fold degenerate hs- 
3d 7 . The spin-orbit state is, as stated above, ls-Zd 7 . Each 
of the manifolds couples to respective hole-transfer states 
which are labeled by Oh irreps of the 3d electron sector 
3 A 2 , 4 Ai, 1 E, 3 Ti, and 3 T 2 of the 3d 8 configuration on 
the cobalt atom. Due to suppression of d 9 occupancy by 
Ud, we restrict the Hilbert space to the one-hole sector 
on the ligands. The states we retain as our basis gives a 
total of 65 states. 

For Co(NH3) 6 + , we consider the following 3<i 6 states: 
The experimentally observed groundstate configuration 
Ai, and the states that couples via spin-orbit coupling 
3 Ti, 3 T2. We would like to observe the suppression of the 
d s configuration, so we include many more hole-transfer 
states for this calculation. For the 3d 7 L states we have 
2 E and 4 Ti. We also include the two- hole states (3d s L 2 ) 
3 A 2 , 1 A 1 , 4 E, 3 T l5 and 3 T 2 , where the hole-degeneracies 
are fully taken into account, giving us 376 states. For 
the case of both the 2+ and 3+ valences, we assume the 



truncation of the basis set is justified by the suppression 
of all other d n configurations by Ud and all other spin 
states by the Hund's rule exchange interaction Jd- 

The basis sets and the quantum weight distribution for 
Co(NH 3 ) 2+ and Co(NH 3 ) 3+ are shown in Fig. (6). The 
next section discusses how parameters were obtained for 
our model. 



C. Determination of model parameters 

Most of our model parameters are extracted from sin- 
gle complex calculations with the fully ab initio code 
SIESTA [13 based on density functional theory (DFT). 
Before discussing the procedure in detail, we briefly want 
to introduce this code. Nevertheless we want to stress 
that our way of determining the parameters is quite gen- 
eral and independent of the code one uses as long as it 
has certain features such as a local atomic basis set and 
the projected density of states (PDOS). 

In running SIESTA we used Troullier-Martins norm- 
conserving pseudo potentials [3l| in the Klcinman- 
Bylander form [33|. For cobalt, we included spin- 
polarization and non-linear core corrections (33J to ac- 
count for a spin-dependent exchange splitting and cor- 
relation effects between core and valence electrons, re- 
spectively. Relativistic effects are included for the core 
electrons in the usual scalar-relativistic approximation 
(mass-velocity and Darwin terms) and by averaging over 
spin-orbit coupling terms, while no spin-orbit coupling 
is included for the 4s and 3d valence electrons. Hence 
spin stays a good quantum number. SIESTA uses a lo- 
cal basis set of pseudo atomic orbitals (PAO) of multiple 
C-type. The first-C orbitals are produced by the method 
by Sankey and Niklewski J34| , while the higher-^ orbitals 
are obtained from the split valence method well known 
from quantum chemistry. Polarization orbitals can also 
be included. We used a double-^ basis set with polariza- 
tion orbitals (DZP), as well as the generalized gradient 
approximation (GGA) in the version by Perdew, Burke 
and Ernzerhof 35] for the exchange-correlation energy 
functional. As for the single complex calculations the 
structures were optimized with the conjugate gradient 
method untill the forces on each atom were below a tol- 
erance of O.OleV/A 

We first discuss the estimation of the following param- 
eters, the d-d Coulomb repulsion Ud, Hund's coupling 
Jd, bare (undressed) d-orbital energy level ej, the ligand 
energy level e/,, the ligand-field (purily due to electro- 
statics) A , and 10D q which contains contributions from 
both electrostatics and hybridization. In order to get a 
physical picture of the single complex problem we show 
the PDOS (DOS projected on atomic 4s, 3 i 2g and 3 e g ) 
of Co(III) in figure |21 as well as Co(II, spin-up) in figure 
Eland Co(II, spin-down) in figure^ All states below the 
Fermi energy tp (set to zero) are occupied. One can im- 
mediately see the 10D q energy gap between t 2g and the 
anti-bonding e* orbitals due to the approximately octra- 



hedral symmetry. The 7r-symmetric t 2g orbitals hardly 
hybridize with the er-symmetric HOMO of NH 3 while 
the cr-symmetric e g orbitals form bonding e g and anti- 
bonding e* molecular orbitals (MOs). For further clarifi- 
cation this is schematically illustrated in figure As for 
Co(II) the non-vanishing total spin (S=3/2) introduces 
spin-dependent MOs due to Hund's coupling, and the 
partly occupied t% g orbitals split probably due to rhom- 
bic distortions. Interestingly, the e g weights are of simi- 
lar magnitude in the bonding e g and the antibonding e* 
MOs. This indicates the near degeneracy of the e g and 
ligand energy level e d (HOMO of NH3). In the case of 
Co (II, spin- up) [see Fig. |3| they are exactly degenerate 
which provides the ligand energy level e^ + = 0.5(e* + e g ) 
immidiately. The t2 g and 10D q energies of Co(II/III) can 
be read off the PDOS, too. In the case of Co(II) we use 
the weighted average of the ti g peaks. The remaining 
parameters are estimated as follows. 

In order to get the ligand energy level e^ + of Co(III) 
we use a simple two-level approximation to describe the 
mixing of the C0-NH3 levels. In matrix form this is 



e g V 
V e L 



(5) 



where V is the hybridization between the Co e g levels 
and the ammine level cl- Setting the two eigenvalues 
of Eq. (0 equal to e g and e* and also constraining the 
coefficients of the eigenvectors by the weights from the 
PDOS in figure |2 we obtain e^ + (as well as V and e g ). 

The Coulomb parameters were obtained from a crude 
mapping assumed between the Kohn-Sham energies ob- 
tained from the GGA calculation and a mean-field calcu- 
lation, where we make use of the single-particle energies 
of the Co t 2g levels 



~df t ^ ~HF 

t t1 g ,o ~ t1 g ,o 



Jd, 



e d + U d (N d -l) + -^(N, 



da 



N d 



-gA (6) 



where e d is the bare energy (without any electron-electron 
interaction), N da are the number of 3d electrons with spin 
a, N d is the total number of d-electrons, and Ao is the 
Madelung (point-charge) contribution to the ligand field. 
We can write down three equations of this kind, one for 
Co(III) (spin-independent) and two for Co(II) (one for 
spin-up and another one for spin-down), and can solve 
for the three unknowns e d , U d , and J d . The required 
parameter A can be otained as follows. 

Since we know the ti g energies from the PDOS we need 
to know the e g energies in order to determine the electro- 
static contribution A = e g — t^g for Co (II) and Co (III). 
These can be extracted in a straightforward way by look- 
ing sharply at figure [S] which summarizes all the involved 
quantities. From the center e = 0.5(e g + e*) of the bond- 
ing e g and anti-bonding e* MOs and the knowledge of 
the ligand energy level e^ we otain e g = e-t- (e — e^). The 
resulting point charge contributions to the total ligand- 
field 10D q are small, A 2 + = l.leV and A^+ = O.leV. 



This is because DFT seems to be biased towards form- 
ing rather extended MOs due to an improper treatment 
of the on-site repulsion among d electrons leading to the 
well-known deficiency of over-binding [7J. It does not 
seem to be surprising that A 2+ > A^ + since due to a 
0.2A shorter Co-N bond length in the case of Co(III) the 
hybridization is even larger than for Co(II). Within DFT 
the lODg of Co(III) is mainly due to hybridization. 

The ability of an electron to transfer from one cobalt 
ion to the other relies on the hybridization between the 
HOMO's on the two bridging ligands, which is estimated 
from the splitting of the HOMO peaks of a two 2 ammine 
calculation. The ammines are in the apex-to-apex orien- 
tation. The splitting of the HOMO corresponds to 2ti nter 
where we define U nter to be the inter complex NH3-NH3 
hopping integral. A similar procedure was applied to ob- 
tain the intra complex NH3-NH3 hopping integral U n tra 
between ammines of the same complex allowing holes to 
delocalize among the ligands. The Co — NH3 hybridiza- 
tion V was determined by computing the expectation val- 
ues <x 2 ~y 2 /z 2 \H\HOMO> for Co(II) and Co(III) us- 
ing the HOMO wavefunction from a single ammine and 
the Hamiltonian matrix of a single complex calculation. 
Contributions to the hybridization from direct overlap 
of non-orthogonal d and HOMO wavefunctions were ne- 
glected. 

Finally, we determine the phonon-coupling constant g 
and the spring constant K as already outlined in section 
ITT1 According to Eq. J3J) g depends on the Slater-Koster 
parameter r\ vdtJ . It is determined by setting V pda = V 
from the DFT method. As for K we use the experimental 
Co-ligand stretching frequencies, 357cm -1 for Co(II) and 
494cm -1 for Co(III) J2J|. 

The parameters are summarized in Table I. Therefore 
we can write down an effective Hamiltonian that includes 
only these frontier orbitals: the set of Co 3d orbitals 
and the NH3 HOMO states. We use the Hamiltonian of 
Eq.(l) as the effective Hamiltonian for each of the sin- 
gle cobalt hexaammine complexes, the 2+ and 3+ nom- 
inal valences. The d-electron parameters e d , U d , J d , the 
ligand-field splitting lODq, etc. all refer to the cobalt 3d 
electrons, while the ligand label L refers to the ammine 
HOMO. We would like to solve for the groundstates of 



Co(NH 3 ) 



2+/3+ 



within this model. 



IV. RESULTS AND DISCUSSION 
A. Quantum weights of single complexes 

The model Hamiltonian we constructed is diagonalized 
on the many-electron basis we have described using the 
set of parameters presented in Table I. We show the re- 
sults in Fig. (6). 

For Co(NH 3 )2+, there is 71.5% in the hs-3d 7 ( 4 Ti) 
state and 28.1% in the corresponding 3d 8 L 2 where an 
electron has transferred from the NH3 HOMO to the e g 
shell of the Co atom. This weight is actually distributed 
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FIG. 2: Projected density-of-states (PDOS) of low-spin 
Co(III) (d , S — 0) for spin-up and spin-down. The Fermi 
energy zf = — 11.82eV is set to zero. 
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FIG. 3: PDOS of high-spin Co(II) (d 7 , S = 3/2). The spin-up 
part is shown. The Fermi energy sf — — 5.77eV is set to zero. 
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FIG. 4: PDOS of high-spin Co(II) (d 7 ,S = 3/2). The spin- 
down part is shown. The Fermi energy 6f = — 8.176 1 !/ is set 
to zero. 
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with ~5% on each of the six coordinated ammine lig- 
ands. This distribution of quantum weight is due to the 
NH3-NH3 hybridization we included on our model. One 
may picture the hole as delocalized or smeared over the 
first solvation shell. The fact that only a small fraction of 
the hole resides on any one ligand has important implica- 
tions for the charge-transfer process and we will discuss 
this later. The remainder of the quantum weight resides 
in the ls-3d 7 ( 2 E), which is small since it mixes into the 
groundstate by the weak spin-orbit coupling. The verti- 
cal energy separation E^E—^Ti) is 1.18 eV (9513 cm -1 ), 
which is in good agreement with ab initio calculations. 

For Co(NH 3 )g + , the vertical energy separation 
E( 3 Ti — 1 Ai) is too small (~0.5 eV) when compared with 



FIG. 5: Ligand-field (LF) and molecular orbital (MO) pic- 
tures. 



experiment (~1.7 eV) for the set of parameters. A pa- 
rameter that we have some freedom with is the ligand- 
field splitting lQDq. The lack of constraint on lODq 
emerges from the model not taking proper account of 
all of the contributions to the ligand-field energy split- 
ting. lQDq is composed of the point-charge ligand-field 
Ao and a contribution due to the hybridization of the 
cobalt and ligand orbitals. While the model does take 
into account the most important hybridization between 
the e g and ammine HOMO, it does not actually form 
bonding/anti-bonding pairs since the 3d electron sector 



TABLE I: Parameters determined for hexaammine cobalt 
complexes from SIESTA calculations. (1) and (2) refer to in- 
plane and axial ligands, respectively. All parameters except 
g-y(eV/A) and K(eV/A 2 ) are in units of eV.) 



Parameters 


Co(II) 


Co(III) 


id 




-33.9 


-33.9 


u d 




4.1 


4.1 


Jd 




-2.0 


-2.0 


g~t 




-3.8 


-5.2 


K 




48.0 


92.0 


t 




0.09 


0.09 


WDq 




1.6 


2.6 


£H 




-11.1 


-14.4 


V(x 2 - 


-y 2 -L) 


1.9 


2.8 


V(z 2 - 


-L){1) 


-1.1 


-1.6 


V(z 2 - 


- L) (2) 


-2.1 


-3.2 



Cu(UI) 




(0.03%) 



FIG. 6: Many-body groundstates of Co(NH 3 )g+ and 
Co(NH 3 )g + . The quantum weights are the squared expan- 
sion coefficients of the groundstate eigenvector | a t \ 2 on the 
restricted Hilbert space. The red arrows indicate states that 
would couple through the bridge Hamiltonian. The blue ar- 
rows show coupling through the spin-orbit term. There are 
65 basis states for Co(NH 3 )g + and 376 for Co(NH 3 )e+. 



of our basis set is represented in the atomic orbital basis, 
i.e. the e g levels are just the atomic orbitals, so in the 
current treatment we cannot accurately predict the hy- 
bridization contribution. Therefore, we treat the ligand- 
field splitting as a free-parameter in our model. While 
the lQDq value seems to work well for Co(NH 3 )g + it does 
not reproduce the correct vertical energy separation be- 
tween low-spin and S — 1 3d 6 , so by increasing it by 1 
eV to 3.55 eV, we recover the vertical energy from ex- 
periment. The effect on the weights is to shift ^5% from 
ls-3d 7 to ls-3d 6 , as would be expected. 

We now need to determine if and how the effects of cor- 
relation between the d electrons manifest in this system. 
For Co(NH 3 )g + , the hole seems nearly fully delocalized 
between the first coordination sphere of ammines and the 
cobalt 3d shell with 51% in the ls-3d 6 ( 1 Ai) configuration 



and 42% in the ls-3d 7 L states. This weight is divided up 
with ~7% per state with a hole localized in one of the six 
NH3 HOMO's. At first glance, the physical picture ap- 
pears to be very much like the screening cloud associated 
with the Kondo effect. In the case of cerium materi- 
als, the hole left behind from an electron hopping from 
the conduction band to the 4/ impurity delocalizes over 
the entire continuum of conduction states, dynamically 
screening the local moment and leading to the formation 
of the singlet bound state. 

We cannot say if the appreciable hole-weight delocal- 
ized over the ligands is due to mixed valence state or 
not. Both the ls-3d 6 ( 1 A 1 ) and Zs-3d 7 ( 2 E) states possess 
filled diamagnetic t^g shells, so the Hund's rule exchange 
contribution to the d-d Coulomb interaction is equal in 
each case. Since we only have two configurations, we may 
shift the Fermi level to correspond with the ls-3d 7 , there- 
fore we do not have to refer to Ud- From this perspec- 
tive, the derealization of the hole between metal 3d and 
the ammine orbital looks simply like the consequence of 
single-particle covalency. This interpretation is also con- 
sistent with the composition of molecular orbitals from 
our ab initio calculations which show nearly equal con- 
tributions from metal and ligand orbitals in the (e g -cr) 
bonding/anti-bonding pair. 

However, there does seem to be an analogy between 
the d 6 , d 7 L, and d 8 L 2 configurations of Co(NH 3 )g + and 
the configurations f°, f 1 ^, and / 2 L 2 of the mixed va- 
lence rare-earth compounds such as Ce02 and metallic 
a-Ce, where mixed valence has been established. In both 
the transition metal and lanthanide cases these multiple 
configurations mix into the groundstate because of their 
near-degeneracy. In both cases there is also large sup- 
pression of the two ligand- hole states due to the Hubbard- 
U Coulomb repulsion of electrons on the localized im- 
purity levels. In the cerium materials the f 1 electron 
and the ligand hole L can form an exciton-like bound- 
state that preserves the symmetry of the diamagnetic f° 
state, i.e. a spin-singlet state with no orbital moment. 
In Co(NH3) 6 , we would have a similar energy lowering 
associated with freezing out the 2-fold orbital degeneracy 
of the e g shell and the spin-degree of freedom, antiferro- 
magnetically coupled to the hole delocalized over the six 
NH 3 HOMO levels. 

Such a state usually forms if the singly-occupied con- 
figuration is slightly lower in energy than the empty 
level. This requirement is not satisfied in Co(NH3)g + 
since the diamagnetic 3d 6 state is stabilized by 1.6 eV 
with respect to the ls-3d 7 configuration. However, the 
hybridization between the e g and ammine HOMO's are 
large enough (~2-3 eV) to mitigate this energy gap. In 
Ce02 the partial occupancy of the 4/-level is observed 
to be ~0.6, consistent with the lattice volumes that are 
somewhat shorter than expected for Ce 3+ compounds. In 



Co(NH 3 ) 6 , if this type of partial occupancy exists (^0.5 
from our model calculation) , experimentalists should ob- 
serve bond lengths that are somewhere intermediate be- 
tween those expected for Co(NH3)g^ 



and Co(NH 3 ) 3+ . 



For Co(NH 3 )g + , as stated above, the values for the 
single-particle Co 3d and ammine HOMO energies are 
such that there is almost complete resonance, and molec- 
ular orbitals are constructed with ^50% of each type of 
orbital. We observe in our calculated result the effect of 
the large Ud in suppressing the 3d 8 weight so the admix- 
ture is now on the order of ~70% of the hole is localized 
in the d shell, while ^30% is delocalized over the am- 
mine ligands. This suppression of 3d 8 configurations is 
expected. The question is what consequence do these ef- 
fects have on the electron transfer matrix element Hda- 



B. Construction of Born-Oppenheimer Potential 
Energy Surfaces and Estimate for Hda 

We performed calculations of the groundstates for the 
single Co(NH 3 )g+ and Co(NH 3 )g + complexes. We would 
now like to use the basis sets of both of these complexes to 
form an encounter complex that represents the reactants 
in the electron-transfer process. 

We connect the single 2+/3+ complexes in the apex- 
to-apex orientation through a bridge composed of am- 
mines from each molecule, which gives us the Hamilto- 
nian 



H = hi + hi - ty^X c LAa c RA<j + h.c.) 



(7) 



Both hi and I12 refer to the identical model Hamilto- 
nian for each single complex given by Eqn.(l). The third 
term describes hole tunneling through the ligand bridge 
composed of the two ammines in close contact. The cre- 
ation operator c LAa . refers to the creation of a cr-spin hole 
on the left NH 3 ligand (belonging to Co(NH 3 ) 3+ ) and 
c RAa ^ S t ne corresponding operator for the right ammine 
(belonging to Co(NH 3 )g + ). The hopping integral t can 
be estimated by noting that It is approximately the en- 
ergy splitting between the bonding/anti-bonding orbitals 
formed from the HOMO of each of the two bridging NH 3 , 
oriented by close contact of the hydrogen atoms of each 
ammine. We must now solve for the groundstate using 
the appropriate basis set. We construct this by taking 
a tensor product of our basis sets for Co(NH 3 )g + and 
Co(NH 3 )g + . The basis set is doubled to account for ex- 
changing Co(NH 3 )g+ and Co(NH 3 )|+ in space. There is 
then another multiplication by two in order to treat the 
Hermitian matrix in order to use diagonalization routines 
suited for real, symmetric matrices [37| |. So the basis set 
is then of rank 



TV = 376 x 65 x 2 x 2 = 97760. 



(8) 



Given that most routines scale as ~ TV 3 , where N is the 
dimensionality of the Hilbcrt space, conventional rou- 
tines, even on very fast machines (~4 GHz micropro- 
cessor) would make single-point calculations very time- 
consuming. Given that we would like to sample the 2D 
groundstate energy surface as a function of the breathing 




FIG. 7: The Born-Oppenheimer potential energy surface for 
the Co(NH3) 6 ' self-exchange electron transfer reaction. 
Only the inner-sphere contribution is shown. The groundstate 
energy was plotted on a 40x40 grid using a Lanczos diagonl- 
ization routine to diagonalize the 97760 2 eigenvalue problem. 
The contour level spacing is chosen to be 0.03 eV/level. xl 
and xr are the reaction coordinates corresponding to the to- 
tally symmetric metal-ligand stretching mode. The transi- 
tion state, defined as the saddle point along x— Xl—xr is at 
x=0A4A. The difference in the Co-N bond length between 
Co(II) and Co(III) is 0.2k. 



coordinates for each complex, xl and xr at the resolu- 
tion needed to accurately determine the transition state, 
we need very efficient matrix diagonalizer. Plotting this 
potential energy surface is a computationally intensive 
problem. 

The solution was to implement a Ldnczos matrix- 
diagonalization routine, which scales as order N [2a| since 
it involves only matrix-vector multiplication for sparse 
matrices. The Ldnczos algorithm converts a large N x N 
matrix into a tridiagonal m x m matrix (m < M). ft 
uses a three-term recurrence relationship which requires 
one trial groundstate vector to start the sequence. The 
eigenvectors tend to become parallel during iteration of 
the procedure. Since we are not interested in them for 
this application, a costly (~ TV 3 ) reorthogonalization step 
does not have to be implemented. 

The plot in Fig. (9) is a two-dimensional contour plot 
representing the Born-Oppenheimer potential energy sur- 
face we have plotted solving the 97760 2 eigenvalue prob- 
lem 1600 times. 

The oval shape of the minima are consistent 
with a larger force constant (narrower potential) for 
Co(NH 3 )g + . The horizontal and vertical distances be- 
tween the minima correspond to ^0.2A which agrees 
well with the experimentally measured difference in the 
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lations. We have essentially reproduced this calculation 
from our simpler model. What the reduction of | H da | 
implies is that the correlation effects we have observed in 
the single complexes of cobalt hexaammine do not seem 
to enhance the ET rate. The reason is that we still rely on 
states that couple through the spin-orbit coupling, most 
noticeably ls-3d 7 ( 2 E) for the Co(NH 3 )j? + complex. The 
quantum weight of this state in the groundstate is not 
increased by the correlation-induced suppresion of 3d 8 
occupancy. In fact, this suppression tends to reduce the 
magnitude of Hda- The arrows in Fig. (2) indicate that 
the 3d 8 states of both Co(NH 3 )^+ and Co(NH 3 )^+ cou- 
ple through the bridge term of H. The weight on these 
states is small because of the large Ud value. 



FIG. 8: Cross-section of PES from Fig. (7) along the reaction 
coordinate d—XL-XR close to the TS. The energy splittings 
(2| Hba I) between the groundstate Eq and first excited state 
Ei at the TS (0A) is about 2 cm -1 . The splitting between 
Ei and E 2 at d=-2.8x!0 _4 A has a similar value. 



metal-ligand distance between the 2+ and 3+ forms of 
the hexaammine cobalt complex. The donor surface 
D for the ET corresponds to the encounter complex 
Co(NH 3 )t + (NH3)-(H 3 N)Co(NH3)5 + , which is given by 



the parabola in the upper left corner. The acceptor 
surface A is that of Co(NH 3 )2+(NH3)-(H 3 N)Co(NH3)^ + 
and is the lower right parabola in Fig. (7) . The transition 
state is the cusp feature at x=0AAA and corresponds to a 
narrowly-avoided crossing between the D and A diabatic 
surfaces. The electronic coupling opens a gap equal to 
2 1 Hda I between D and A at the crossing point, where 
Hda is the electronic coupling factor. 

The general method of extracting Hda is as follows: 
We find the TS by plotting the energy along the diag- 
onal x in Fig. (7) with high resolution. The TS is the 
saddlepoint at x=0A4A. We then read off the difference 
between the first excited and ground states. Hda is half 
of this difference. 

On the other hand, one can exploit that Hda is also a 
groundstate property. Instead of calculating the exictcd 
state, one can obtain it from calculating the energy dif- 
ference between the groundstate adiabatic surface eval- 
uated at the TS and the avoided crossing between the 
diabatic D and A surfaces that results when the bridge 
hopping integral t is turned off. 

In Fig. (10) we plot a cross-section of PES from Fig. (9) 
along the reaction coordinate d — x^-x^ in the region 
of the TS. The energy splittings (2| Hda |) between the 
groundstate Eq and first excited state E\ at the TS (0A) 
is ~2 cm -1 . The splitting between E\ and Ei at d=- 
2.8xlCP 4 A has a similar value. The splittings between 
Eq and Ei and also E\ and E^ originate from lifting of 
the degeneracy of the components of the 4 Ti irreducible 
representation of Co(II) by the spin-orbit coupling. 

Our result (<~2 cm -1 ) is smaller than Newton's calcu- 
lation of Hda, which gave 9.5 cm -1 from ab initio calcu- 



V. CONCLUSIONS 

We have introduced a variant of the Anderson Im- 
purity Model that is appropriate for describing the va- 
lence electronic structure of TM molecules. We have 
coupled this model to a metal-ligand stretching mode. 
We have applied this model to investigate the effect of 
strong correlations among the cobalt 3d electrons on the 
self-exchange electron transfer rate in the cobalt hexam- 
mines Co(NH 3 )g +/3+ . 

We have constrained the model parameters by using 
the results of ab initio density functional theory cal- 
culations and experimental stretching frequencies. The 
single-particle hybridizations VdL were extracted directly 
and the Coulomb parameters (ed,Ud,Jd) were obtained 
from mapping to a mean-field solution of our model. Us- 
ing these parameters, we found the groundstates of both 
2+ and 3+ nominal valences. For the Co(NH 3 )g + state 
we find strong configurational admixture between the 3d 6 
and 3d 7 L states. The hole L is fully delocalized between 
the cobalt atom and the first coordination sphere of am- 
mines. While this resembles the screening cloud and par- 
tial occupancy of the impurity level seen in mixed valcnt 
compounds like Ce02, the large 3d 7 L weight could simply 
be a purely single-particle effect since the derealization 
is consistent with the equal metal and ligand orbital com- 
position of molecular orbitals near the Fermi level from 
the DFT calculation. For Co(NH3)g + there is an obvious 
signature of the effect of large Ud, with a suppression of 
~20% of the quantum weight on 3d 8 L 2 states. 

We couple the complexes together by a bridging Hamil- 
tonian to form the encounter complex and generate a 
basis set by using the tensor product of the single hex- 
aammine basis states. This gives a very large eigenvalue 
problem (N ~10 5 ) that is solved using Ldnczos meth- 
ods. The groundstate energy surface is plotted, resolving 
the two minima of the electron transfer process. We are 
able to extract the electronic coupling Hda as the energy 
lowering at the transition state, relative to the diabatic 
surfaces. What we found is that configurational admix- 
ture does not enhance Hda and we get a result ^4 cm -1 , 
that is within a factor of 5 of Newton's original calcula- 
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tion. 

It appears that the disparity of two orders of mag- 
nitude in the self-exchange electron transfer rate of 
Co(NH3) 6 cannot be explained on the basis of the 

electronic structure of a single complex. One aspect that 
has been completely ignored in this calculation are sol- 
vation effects beyond the first coordination shell. While 
one would expect, in the gas phase, that the single, sym- 
metric breathing coordinate would be an excellent ap- 
proximation to the breathing coordinate, solvent effects 
could significantly change the heigth of the barrier by 
adding anharmonicity to the potential wells. This could 
contribute to the adiabatic nature of the reaction. 

Another possible source of enhancement for the ET 
rate could be through dynamical effects. The dynami- 
cal prefactor depends on the details of how the system 
crosses and recrosses the thermal barrier. The presence 
of the low-lying spin-orbit states that are separated from 
the groundstate by very small (~1 cm" 1 ) energies could 
give the vibronic wavepacket multiple channels for reach- 
ing the other minimum. This is a possibility that can be 
explored with our model, which provides fast "on-the- 
fly" sampling of the correct electronic structure. Excited 
state calculations of the bimetallic TM active site involv- 
ing open-shell states, coupled by the TM atom spin-orbit 
interaction, would be very difficult to do with state-of- 
the-art quantum chemistry codes. 

Finally, it is possible that the relevant Born- 



Oppenheimer surfaces at the transition state are differ- 
ent from those at equilibrium, which we have explored 
elsewhere within DFT theory [l^J. In particular, an in- 
termediate spin surface appears to be thermally accessi- 
ble and favorable at the transition state. In that paper 
we actually found a rate two orders of magnitude faster 
than experiment. It is interesting to note that within 
our calculation and Newton's UMP2 calculationi 18] the 
electron correlation effects actually reduce the transition 
rate. Within our variational state treatment the reason 
is clear: Coulombic blocking effects eliminate a class of 
determinants which can facilitate charge transfer. It will 
be interesting to explore the effects of correlations upon 
the intermediate spin manifold of states. 

The method proposed in this paper could be applied 
to other TM systems such as the active sites in metallo- 
proteins. More elaborate ways of coupling the electronic 
sector to the protein structure will have to be devised 
but the efficiency and simplicity of our method might be 
very useful in identifying transition state configurations. 
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